Just in case we don’t have our data saved from the last round:
library(dismo)
## Loading required package: raster
## Loading required package: sp
library(rgeos)
## rgeos version: 0.5-3, (SVN revision 634)
## GEOS runtime version: 3.7.2-CAPI-1.11.2
## Linking to sp version: 1.4-1
## Polygon checking: TRUE
library(rgbif)
library(leaflet)
ibm <- occ_search(scientificName = "Iberolacerta monticola", limit = 1500)
keep.cols <- c("species", "decimalLatitude", "decimalLongitude")
ibm <- ibm$data[,keep.cols]
ibm <- as.data.frame(unique(ibm))
ibm <- ibm[complete.cases(ibm),]
colnames(ibm) <- c("species", "lat", "lon")
head(ibm)
## species lat lon
## 1 Iberolacerta monticola 43.45172 -6.699292
## 2 Iberolacerta monticola 43.43051 -7.188726
## 3 Iberolacerta monticola 43.46805 -6.703289
## 4 Iberolacerta monticola 43.51037 -6.504768
## 5 Iberolacerta monticola 43.52801 -6.556106
## 6 Iberolacerta monticola 43.55888 -6.728039
all.worldclim <- raster::getData("worldclim", res = 10, var = "bio")
spain.worldclim <- crop(all.worldclim, extent(-10, 4, 35, 45))
m <- leaflet(ibm) %>%
addProviderTiles(provider = "Stamen.TerrainBackground") %>%
addRasterImage(spain.worldclim[[1]]) %>%
addCircleMarkers(lng = ~lon, lat = ~lat, radius = 5, stroke = FALSE,
fillOpacity = 1)
## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO"): Discarded
## ellps WGS 84 in CRS definition: +proj=merc +a=6378137 +b=6378137 +lat_ts=0
## +lon_0=0 +x_0=0 +y_0=0 +k=1 +units=m +nadgrids=@null +wktext +no_defs
## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO"): Discarded
## datum WGS_1984 in CRS definition
## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO"): Discarded
## ellps WGS 84 in CRS definition: +proj=merc +a=6378137 +b=6378137 +lat_ts=0
## +lon_0=0 +x_0=0 +y_0=0 +k=1 +units=m +nadgrids=@null +wktext +no_defs
## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO"): Discarded
## datum WGS_1984 in CRS definition
m
Oops! We’ve got one point down there on the west coast of Africa. We could explore why that is, or we could just eliminate it. For the sake of brevity, we’ll do the latter.
ibm <- ibm %>%
dplyr::filter(lat > 30)
m <- leaflet(ibm) %>%
addProviderTiles(provider = "Stamen.TerrainBackground") %>%
addRasterImage(spain.worldclim[[1]]) %>%
addCircleMarkers(lng = ~lon, lat = ~lat, radius = 5, stroke = FALSE,
fillOpacity = 1)
## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO"): Discarded
## ellps WGS 84 in CRS definition: +proj=merc +a=6378137 +b=6378137 +lat_ts=0
## +lon_0=0 +x_0=0 +y_0=0 +k=1 +units=m +nadgrids=@null +wktext +no_defs
## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO"): Discarded
## datum WGS_1984 in CRS definition
## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO"): Discarded
## ellps WGS 84 in CRS definition: +proj=merc +a=6378137 +b=6378137 +lat_ts=0
## +lon_0=0 +x_0=0 +y_0=0 +k=1 +units=m +nadgrids=@null +wktext +no_defs
## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO"): Discarded
## datum WGS_1984 in CRS definition
m
Nice!
Now we’ll load ENMTools. If you’ve already installed it and run install.extras(), you can run it just by typing:
library(ENMTools)
If you haven’t done either of those things, see the “Installing ENMTools” tutorial for details.
ENMTools automates a lot of the fiddly parts of building SDMs, but in order to do that it requires us to get our data into a specific format.
monticola <- enmtools.species()
monticola
##
##
## Range raster not defined.
##
## Presence points not defined.
##
## Background points not defined.
##
## Species name not defined.
We now have an empty enmtools.species object. It’s basically just a list where we can store our data. This makes sure that our data is in a predictable format for ENMTools’ modeling functions. Let’s go ahead and add data to it now.
monticola$presence.points <- ibm[,c("lon", "lat")]
monticola$species.name <- "Iberolacerta monticola"
monticola
##
##
## Range raster not defined.
##
## Presence points (first ten only):
##
## | lon| lat|
## |---------:|--------:|
## | -6.699292| 43.45172|
## | -7.188726| 43.43051|
## | -6.703289| 43.46805|
## | -6.504768| 43.51037|
## | -6.556106| 43.52801|
## | -6.728039| 43.55888|
## | -6.823801| 43.39811|
## | -6.732560| 43.54083|
## | -6.650406| 43.40437|
## | -6.541960| 43.45676|
##
##
## Background points not defined.
##
## Species name: Iberolacerta monticola
We can define a range for our species using any raster, or we can build one using the background.raster.buffer function of ENMTools
monticola$range <- background.raster.buffer(monticola$presence.points, 50000, mask = spain.worldclim)
plot(monticola)
This range will be used for sampling background points in ENM construction, but we could also provide points manually if we wanted to using monticola$background.points <- something.
It’s always a good idea to run your species through the “check.species” function to make sure it’s formatted correctly
monticola <- check.species(monticola)
Not only does this tell us if we’ve messed something up, it reformats column names etc. so that they’re in the appropriate format for the remaining ENMTools features. The fact that everything is stored in a reliable place in these enmtools.species objects allows us to streamline a lot of operations. For instance, instead of manually building a leaflet plot we can now just use the enmtools.species object to build one automatically.
interactive.plot.enmtools.species(monticola)
## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO"): Discarded
## ellps WGS 84 in CRS definition: +proj=merc +a=6378137 +b=6378137 +lat_ts=0
## +lon_0=0 +x_0=0 +y_0=0 +k=1 +units=m +nadgrids=@null +wktext +no_defs
## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO"): Discarded
## datum WGS_1984 in CRS definition
## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO"): Discarded
## ellps WGS 84 in CRS definition: +proj=merc +a=6378137 +b=6378137 +lat_ts=0
## +lon_0=0 +x_0=0 +y_0=0 +k=1 +units=m +nadgrids=@null +wktext +no_defs
## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO"): Discarded
## datum WGS_1984 in CRS definition
Okay, we’ve got a species. Let’s build a model! We’ll start with the two most popular climate envelope models.
monticola.bc <- enmtools.bc(species = monticola, env = spain.worldclim, test.prop = 0.3, nback = 400)
monticola.bc
##
##
## Data table (top ten lines):
##
## | | Longitude| Latitude|
## |:--|---------:|--------:|
## |1 | -6.699292| 43.45172|
## |2 | -7.188726| 43.43051|
## |5 | -6.556106| 43.52801|
## |7 | -6.823801| 43.39811|
## |8 | -6.732560| 43.54083|
## |9 | -6.650406| 43.40437|
## |11 | -6.910225| 43.29794|
## |12 | -6.659564| 43.48790|
## |13 | -6.817537| 43.28005|
## |14 | -6.831939| 43.38955|
##
##
## Model: class : Bioclim
##
## variables: bio1 bio2 bio3 bio4 bio5 bio6 bio7 bio8 bio9 bio10 bio11 bio12 bio13 bio14 bio15 bio16 bio17 bio18 bio19
##
##
## presence points: 294
## bio1 bio2 bio3 bio4 bio5 bio6 bio7 bio8 bio9 bio10 bio11 bio12 bio13 bio14
## 1 121 88 41 4237 239 28 211 97 175 178 70 920 121 36
## 2 128 84 41 4080 239 38 201 83 180 182 78 959 124 33
## 3 134 85 42 3987 243 45 198 113 185 187 86 860 116 35
## 4 121 88 41 4237 239 28 211 97 175 178 70 920 121 36
## 5 137 83 42 3957 244 48 196 115 187 190 89 852 114 33
## 6 117 90 42 4292 237 23 214 93 172 174 65 926 123 38
## 7 117 91 41 4472 240 20 220 67 175 176 63 942 122 34
## 8 117 90 42 4292 237 23 214 93 172 174 65 926 123 38
## 9 105 93 41 4614 233 7 226 54 165 166 49 983 127 38
## 10 121 88 41 4237 239 28 211 97 175 178 70 920 121 36
## bio15 bio16 bio17 bio18 bio19
## 1 31 320 135 148 279
## 2 34 345 126 142 310
## 3 30 299 130 145 248
## 4 31 320 135 148 279
## 5 31 298 124 139 253
## 6 30 320 140 153 274
## 7 33 334 129 142 300
## 8 30 320 140 153 274
## 9 32 345 141 152 309
## 10 31 320 135 148 279
## (... ... ...)
##
##
##
## Model fit (training data): class : ModelEvaluation
## n presences : 294
## n absences : 400
## AUC : 0.8389031
## cor : 0.4607588
## max TPR+TNR at : 0.05092041
##
##
## Environment space model fit (training data): class : ModelEvaluation
## n presences : 294
## n absences : 10000
## AUC : 0.9681876
## cor : 0.7088733
## max TPR+TNR at : 0.0169068
##
##
## Proportion of data wittheld for model fitting: 0.3
##
## Model fit (test data): class : ModelEvaluation
## n presences : 127
## n absences : 400
## AUC : 0.8252657
## cor : 0.4288314
## max TPR+TNR at : 0.05092041
##
##
## Environment space model fit (test data): class : ModelEvaluation
## n presences : 127
## n absences : 10000
## AUC : 0.955861
## cor : 0.6948817
## max TPR+TNR at : 0.02030816
##
##
## Suitability:
## class : RasterLayer
## dimensions : 60, 84, 5040 (nrow, ncol, ncell)
## resolution : 0.1666667, 0.1666667 (x, y)
## extent : -10, 4, 35, 45 (xmin, xmax, ymin, ymax)
## crs : +proj=longlat +datum=WGS84 +no_defs
## source : memory
## names : layer
## values : 0.003401361, 0.6564626 (min, max)
##
##
##
## Notes:
## NULL
monticola.bc$response.plots[3]
## $bio3
monticola.bc$response.plots[8]
## $bio8
visualize.enm(monticola.bc, spain.worldclim, layers = c("bio3", "bio8"))
## $background.plot
## [1] NA
##
## $suit.plot
monticola.dm <- enmtools.dm(species = monticola, env = spain.worldclim, test.prop = 0.3, nback = 400)
monticola.dm
##
##
## Data table (top ten lines):
##
## | | Longitude| Latitude|
## |:--|---------:|--------:|
## |1 | -6.699292| 43.45172|
## |2 | -7.188726| 43.43051|
## |3 | -6.703289| 43.46805|
## |4 | -6.504768| 43.51037|
## |5 | -6.556106| 43.52801|
## |6 | -6.728039| 43.55888|
## |7 | -6.823801| 43.39811|
## |8 | -6.732560| 43.54083|
## |13 | -6.817537| 43.28005|
## |14 | -6.831939| 43.38955|
##
##
## Model: class : Domain
##
## variables: bio1 bio2 bio3 bio4 bio5 bio6 bio7 bio8 bio9 bio10 bio11 bio12 bio13 bio14 bio15 bio16 bio17 bio18 bio19
##
##
## presence points: 294
## bio1 bio2 bio3 bio4 bio5 bio6 bio7 bio8 bio9 bio10 bio11 bio12 bio13 bio14
## 1 121 88 41 4237 239 28 211 97 175 178 70 920 121 36
## 2 128 84 41 4080 239 38 201 83 180 182 78 959 124 33
## 3 121 88 41 4237 239 28 211 97 175 178 70 920 121 36
## 4 134 85 42 3987 243 45 198 113 185 187 86 860 116 35
## 5 134 85 42 3987 243 45 198 113 185 187 86 860 116 35
## 6 137 83 42 3957 244 48 196 115 187 190 89 852 114 33
## 7 121 88 41 4237 239 28 211 97 175 178 70 920 121 36
## 8 137 83 42 3957 244 48 196 115 187 190 89 852 114 33
## 9 105 93 41 4614 233 7 226 54 165 166 49 983 127 38
## 10 121 88 41 4237 239 28 211 97 175 178 70 920 121 36
## bio15 bio16 bio17 bio18 bio19
## 1 31 320 135 148 279
## 2 34 345 126 142 310
## 3 31 320 135 148 279
## 4 30 299 130 145 248
## 5 30 299 130 145 248
## 6 31 298 124 139 253
## 7 31 320 135 148 279
## 8 31 298 124 139 253
## 9 32 345 141 152 309
## 10 31 320 135 148 279
## (... ... ...)
##
##
##
## Model fit (training data): class : ModelEvaluation
## n presences : 294
## n absences : 400
## AUC : 0.7351148
## cor : 0.3874604
## max TPR+TNR at : 0.6216023
##
##
## Environment space model fit (training data): class : ModelEvaluation
## n presences : 294
## n absences : 10000
## AUC : 0.8966745
## cor : 0.4572874
## max TPR+TNR at : 0.5652725
##
##
## Proportion of data wittheld for model fitting: 0.3
##
## Model fit (test data): class : ModelEvaluation
## n presences : 127
## n absences : 400
## AUC : 0.7647441
## cor : 0.390142
## max TPR+TNR at : 0.6216023
##
##
## Environment space model fit (test data): class : ModelEvaluation
## n presences : 127
## n absences : 10000
## AUC : 0.9254134
## cor : 0.3550469
## max TPR+TNR at : 0.5534922
##
##
## Suitability:
## class : RasterLayer
## dimensions : 60, 84, 5040 (nrow, ncol, ncell)
## resolution : 0.1666667, 0.1666667 (x, y)
## extent : -10, 4, 35, 45 (xmin, xmax, ymin, ymax)
## crs : +proj=longlat +datum=WGS84 +no_defs
## source : memory
## names : layer
## values : 0.3631167, 0.7752813 (min, max)
##
##
##
## Notes:
## NULL
monticola.dm$response.plots[3]
## $bio3
monticola.dm$response.plots[8]
## $bio8
visualize.enm(monticola.dm, spain.worldclim, layers = c("bio3", "bio8"))
## $background.plot
## [1] NA
##
## $suit.plot
monticola.glm <- enmtools.glm(species = monticola, env = spain.worldclim)
## Warning in predict.lm(object, newdata, se.fit, scale = 1, type = if (type == :
## prediction from a rank-deficient fit may be misleading
## Warning in predict.lm(object, newdata, se.fit, scale = 1, type = if (type == :
## prediction from a rank-deficient fit may be misleading
## Warning in predict.lm(object, newdata, se.fit, scale = 1, type = if (type == :
## prediction from a rank-deficient fit may be misleading
## Warning in predict.lm(object, newdata, se.fit, scale = 1, type = if (type == :
## prediction from a rank-deficient fit may be misleading
## Warning in predict.lm(object, newdata, se.fit, scale = 1, type = if (type == :
## prediction from a rank-deficient fit may be misleading
## Warning in predict.lm(object, newdata, se.fit, scale = 1, type = if (type == :
## prediction from a rank-deficient fit may be misleading
## Warning in predict.lm(object, newdata, se.fit, scale = 1, type = if (type == :
## prediction from a rank-deficient fit may be misleading
## Warning in predict.lm(object, newdata, se.fit, scale = 1, type = if (type == :
## prediction from a rank-deficient fit may be misleading
## Warning in predict.lm(object, newdata, se.fit, scale = 1, type = if (type == :
## prediction from a rank-deficient fit may be misleading
## Warning in predict.lm(object, newdata, se.fit, scale = 1, type = if (type == :
## prediction from a rank-deficient fit may be misleading
## Warning in predict.lm(object, newdata, se.fit, scale = 1, type = if (type == :
## prediction from a rank-deficient fit may be misleading
## Warning in predict.lm(object, newdata, se.fit, scale = 1, type = if (type == :
## prediction from a rank-deficient fit may be misleading
## Warning in predict.lm(object, newdata, se.fit, scale = 1, type = if (type == :
## prediction from a rank-deficient fit may be misleading
## Warning in predict.lm(object, newdata, se.fit, scale = 1, type = if (type == :
## prediction from a rank-deficient fit may be misleading
## Warning in predict.lm(object, newdata, se.fit, scale = 1, type = if (type == :
## prediction from a rank-deficient fit may be misleading
## Warning in predict.lm(object, newdata, se.fit, scale = 1, type = if (type == :
## prediction from a rank-deficient fit may be misleading
## Warning in predict.lm(object, newdata, se.fit, scale = 1, type = if (type == :
## prediction from a rank-deficient fit may be misleading
## Warning in predict.lm(object, newdata, se.fit, scale = 1, type = if (type == :
## prediction from a rank-deficient fit may be misleading
## Warning in predict.lm(object, newdata, se.fit, scale = 1, type = if (type == :
## prediction from a rank-deficient fit may be misleading
## Warning in predict.lm(object, newdata, se.fit, scale = 1, type = if (type == :
## prediction from a rank-deficient fit may be misleading
## Warning in predict.lm(object, newdata, se.fit, scale = 1, type = if (type == :
## prediction from a rank-deficient fit may be misleading
## Warning in predict.lm(object, newdata, se.fit, scale = 1, type = if (type == :
## prediction from a rank-deficient fit may be misleading
## Warning in predict.lm(object, newdata, se.fit, scale = 1, type = if (type == :
## prediction from a rank-deficient fit may be misleading
## Warning in predict.lm(object, newdata, se.fit, scale = 1, type = if (type == :
## prediction from a rank-deficient fit may be misleading
## Warning in predict.lm(object, newdata, se.fit, scale = 1, type = if (type == :
## prediction from a rank-deficient fit may be misleading
## Warning in predict.lm(object, newdata, se.fit, scale = 1, type = if (type == :
## prediction from a rank-deficient fit may be misleading
## Warning in predict.lm(object, newdata, se.fit, scale = 1, type = if (type == :
## prediction from a rank-deficient fit may be misleading
## Warning in predict.lm(object, newdata, se.fit, scale = 1, type = if (type == :
## prediction from a rank-deficient fit may be misleading
## Warning in predict.lm(object, newdata, se.fit, scale = 1, type = if (type == :
## prediction from a rank-deficient fit may be misleading
## Warning in predict.lm(object, newdata, se.fit, scale = 1, type = if (type == :
## prediction from a rank-deficient fit may be misleading
## Warning in predict.lm(object, newdata, se.fit, scale = 1, type = if (type == :
## prediction from a rank-deficient fit may be misleading
warnings()
We’ve got a warning here, which is due to the same reason as before - too many predictors and not enough data.
Let’s look at the model first and then we’ll work on that.
monticola.glm
##
##
## Formula: presence ~ bio1 + bio2 + bio3 + bio4 + bio5 + bio6 + bio7 + bio8 +
## bio9 + bio10 + bio11 + bio12 + bio13 + bio14 + bio15 + bio16 +
## bio17 + bio18 + bio19
## <environment: 0x7fc6ca8ef200>
##
##
## Data table (top ten lines):
##
## | Longitude| Latitude| bio1| bio2| bio3| bio4| bio5| bio6| bio7| bio8| bio9| bio10| bio11| bio12| bio13| bio14| bio15| bio16| bio17| bio18| bio19| presence|
## |---------:|--------:|----:|----:|----:|----:|----:|----:|----:|----:|----:|-----:|-----:|-----:|-----:|-----:|-----:|-----:|-----:|-----:|-----:|--------:|
## | -6.699292| 43.45172| 121| 88| 41| 4237| 239| 28| 211| 97| 175| 178| 70| 920| 121| 36| 31| 320| 135| 148| 279| 1|
## | -7.188726| 43.43051| 128| 84| 41| 4080| 239| 38| 201| 83| 180| 182| 78| 959| 124| 33| 34| 345| 126| 142| 310| 1|
## | -6.703289| 43.46805| 121| 88| 41| 4237| 239| 28| 211| 97| 175| 178| 70| 920| 121| 36| 31| 320| 135| 148| 279| 1|
## | -6.504768| 43.51037| 134| 85| 42| 3987| 243| 45| 198| 113| 185| 187| 86| 860| 116| 35| 30| 299| 130| 145| 248| 1|
## | -6.556106| 43.52801| 134| 85| 42| 3987| 243| 45| 198| 113| 185| 187| 86| 860| 116| 35| 30| 299| 130| 145| 248| 1|
## | -6.728039| 43.55888| 137| 83| 42| 3957| 244| 48| 196| 115| 187| 190| 89| 852| 114| 33| 31| 298| 124| 139| 253| 1|
## | -6.823801| 43.39811| 121| 88| 41| 4237| 239| 28| 211| 97| 175| 178| 70| 920| 121| 36| 31| 320| 135| 148| 279| 1|
## | -6.732560| 43.54083| 137| 83| 42| 3957| 244| 48| 196| 115| 187| 190| 89| 852| 114| 33| 31| 298| 124| 139| 253| 1|
## | -6.650406| 43.40437| 117| 90| 42| 4292| 237| 23| 214| 93| 172| 174| 65| 926| 123| 38| 30| 320| 140| 153| 274| 1|
## | -6.541960| 43.45676| 117| 90| 42| 4292| 237| 23| 214| 93| 172| 174| 65| 926| 123| 38| 30| 320| 140| 153| 274| 1|
##
##
## Model:
## Call:
## glm(formula = f, family = "binomial", data = analysis.df[, -c(1,
## 2)], weights = weights)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -1.5150 -0.6744 -0.2519 0.6263 3.6446
##
## Coefficients: (1 not defined because of singularities)
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -4.476680 13.558908 -0.330 0.741275
## bio1 -0.743396 0.205341 -3.620 0.000294 ***
## bio2 0.252844 0.155605 1.625 0.104182
## bio3 0.462477 0.270438 1.710 0.087246 .
## bio4 0.012248 0.005495 2.229 0.025835 *
## bio5 -0.315935 0.097055 -3.255 0.001133 **
## bio6 0.108148 0.105166 1.028 0.303785
## bio7 NA NA NA NA
## bio8 -0.035162 0.008867 -3.965 7.33e-05 ***
## bio9 0.075749 0.200978 0.377 0.706245
## bio10 0.121069 0.213364 0.567 0.570422
## bio11 0.715685 0.211397 3.386 0.000710 ***
## bio12 -0.014210 0.014713 -0.966 0.334136
## bio13 0.029987 0.041034 0.731 0.464911
## bio14 -0.357382 0.147380 -2.425 0.015313 *
## bio15 -0.274822 0.079348 -3.463 0.000533 ***
## bio16 0.082386 0.030383 2.712 0.006696 **
## bio17 0.051112 0.064050 0.798 0.424869
## bio18 0.006822 0.034196 0.200 0.841861
## bio19 -0.038979 0.027766 -1.404 0.160368
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 1167.26 on 1420 degrees of freedom
## Residual deviance: 807.49 on 1402 degrees of freedom
## AIC: 429.02
##
## Number of Fisher Scoring iterations: 10
##
##
##
## Model fit (training data): class : ModelEvaluation
## n presences : 421
## n absences : 1000
## AUC : 0.8457304
## cor : 0.4492968
## max TPR+TNR at : 0.2591298
##
##
## Environment space model fit (training data): class : ModelEvaluation
## n presences : 421
## n absences : 10000
## AUC : 0.555324
## cor : 0.096314
## max TPR+TNR at : 0.117819
##
##
## Proportion of data wittheld for model fitting: 0
##
## Model fit (test data): [1] NA
##
##
## Environment space model fit (test data): [1] NA
##
##
## Suitability:
## class : RasterLayer
## dimensions : 60, 84, 5040 (nrow, ncol, ncell)
## resolution : 0.1666667, 0.1666667 (x, y)
## extent : -10, 4, 35, 45 (xmin, xmax, ymin, ymax)
## crs : +proj=longlat +datum=WGS84 +no_defs
## source : memory
## names : layer
## values : 3.509868e-09, 0.9999487 (min, max)
##
##
##
## Notes:
names(monticola.glm)
## [1] "species.name" "formula"
## [3] "analysis.df" "test.data"
## [5] "test.prop" "model"
## [7] "training.evaluation" "test.evaluation"
## [9] "env.training.evaluation" "env.test.evaluation"
## [11] "rts.test" "suitability"
## [13] "clamping.strength" "call"
## [15] "notes" "response.plots"
Look at all of that stuff!
monticola.glm$response.plots
## $bio1
##
## $bio2
##
## $bio3
##
## $bio4
##
## $bio5
##
## $bio6
##
## $bio7
##
## $bio8
##
## $bio9
##
## $bio10
##
## $bio11
##
## $bio12
##
## $bio13
##
## $bio14
##
## $bio15
##
## $bio16
##
## $bio17
##
## $bio18
##
## $bio19
We have interactive plots for models too!
interactive.plot.enmtools.model(monticola.glm)
## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO"): Discarded
## ellps WGS 84 in CRS definition: +proj=merc +a=6378137 +b=6378137 +lat_ts=0
## +lon_0=0 +x_0=0 +y_0=0 +k=1 +units=m +nadgrids=@null +wktext +no_defs
## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO"): Discarded
## datum WGS_1984 in CRS definition
## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO"): Discarded
## ellps WGS 84 in CRS definition: +proj=merc +a=6378137 +b=6378137 +lat_ts=0
## +lon_0=0 +x_0=0 +y_0=0 +k=1 +units=m +nadgrids=@null +wktext +no_defs
## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO"): Discarded
## datum WGS_1984 in CRS definition
Let’s try fitting a simpler model
First let’s look at the correlations between variables
raster.cor.matrix(spain.worldclim)
## bio1 bio2 bio3 bio4 bio5 bio6
## bio1 1.000000000 -0.05192243 -0.02536989 -0.06482844 0.7309560 0.89386921
## bio2 -0.051922434 1.00000000 0.25804128 0.68346008 0.5595536 -0.41646754
## bio3 -0.025369890 0.25804128 1.00000000 -0.50487533 -0.1563297 0.09040494
## bio4 -0.064828440 0.68346008 -0.50487533 1.00000000 0.5735377 -0.47341350
## bio5 0.730955951 0.55955362 -0.15632969 0.57353766 1.0000000 0.37935531
## bio6 0.893869213 -0.41646754 0.09040494 -0.47341350 0.3793553 1.00000000
## bio7 -0.031962836 0.88012003 -0.22456546 0.94036051 0.6494830 -0.45715397
## bio8 0.404382488 -0.11713369 -0.09355248 -0.03873301 0.2024989 0.33552051
## bio9 0.683287223 0.20191286 0.01963002 0.11583860 0.6773847 0.54402099
## bio10 0.922563848 0.20792847 -0.22821990 0.32228650 0.9151040 0.66801013
## bio11 0.946033848 -0.27710014 0.11951392 -0.37865009 0.4946439 0.98216169
## bio12 -0.392848782 -0.45597177 0.25444757 -0.57123600 -0.6674771 -0.06869193
## bio13 -0.120028527 -0.51403698 0.25059275 -0.64521027 -0.4924741 0.20166588
## bio14 -0.758765734 -0.25110222 0.10775308 -0.24002532 -0.8123826 -0.56244009
## bio15 0.773143662 -0.07304360 -0.01876198 -0.12025929 0.5822825 0.74503621
## bio16 -0.087315073 -0.49267982 0.26595406 -0.63962927 -0.4507268 0.23058503
## bio17 -0.739725188 -0.26496834 0.12139640 -0.26115476 -0.8098741 -0.53463352
## bio18 -0.713437527 -0.33085626 0.08080540 -0.28904180 -0.8266616 -0.49637460
## bio19 -0.002591276 -0.40575986 0.29614464 -0.58921379 -0.3270008 0.28232294
## bio7 bio8 bio9 bio10 bio11 bio12
## bio1 -0.03196284 0.40438249 0.68328722 0.9225638 0.94603385 -0.39284878
## bio2 0.88012003 -0.11713369 0.20191286 0.2079285 -0.27710014 -0.45597177
## bio3 -0.22456546 -0.09355248 0.01963002 -0.2282199 0.11951392 0.25444757
## bio4 0.94036051 -0.03873301 0.11583860 0.3222865 -0.37865009 -0.57123600
## bio5 0.64948302 0.20249891 0.67738465 0.9151040 0.49464394 -0.66747713
## bio6 -0.45715397 0.33552051 0.54402099 0.6680101 0.98216169 -0.06869193
## bio7 1.00000000 -0.08108260 0.20404940 0.3306597 -0.33167441 -0.58515377
## bio8 -0.08108260 1.00000000 -0.04153511 0.3626530 0.38318128 -0.42603645
## bio9 0.20404940 -0.04153511 1.00000000 0.7008436 0.60235162 -0.31651363
## bio10 0.33065966 0.36265299 0.70084359 1.0000000 0.75367532 -0.59770879
## bio11 -0.33167441 0.38318128 0.60235162 0.7536753 1.00000000 -0.18778961
## bio12 -0.58515377 -0.42603645 -0.31651363 -0.5977088 -0.18778961 1.00000000
## bio13 -0.63911541 -0.33110726 -0.12411293 -0.3619468 0.09551341 0.92211902
## bio14 -0.31867777 -0.22120096 -0.66922270 -0.8260330 -0.64267656 0.62160840
## bio15 -0.05256181 0.06324795 0.60151828 0.7025760 0.77296096 -0.20402064
## bio16 -0.62275227 -0.37105106 -0.08385793 -0.3281109 0.12411688 0.91197261
## bio17 -0.33911807 -0.22760761 -0.63176811 -0.8165397 -0.61836117 0.66182192
## bio18 -0.38669626 -0.13863713 -0.68989269 -0.8016505 -0.58312240 0.64972769
## bio19 -0.54634045 -0.44956955 0.05633320 -0.2288494 0.18398498 0.85972346
## bio13 bio14 bio15 bio16 bio17 bio18
## bio1 -0.12002853 -0.7587657 0.77314366 -0.08731507 -0.7397252 -0.7134375
## bio2 -0.51403698 -0.2511022 -0.07304360 -0.49267982 -0.2649683 -0.3308563
## bio3 0.25059275 0.1077531 -0.01876198 0.26595406 0.1213964 0.0808054
## bio4 -0.64521027 -0.2400253 -0.12025929 -0.63962927 -0.2611548 -0.2890418
## bio5 -0.49247406 -0.8123826 0.58228255 -0.45072685 -0.8098741 -0.8266616
## bio6 0.20166588 -0.5624401 0.74503621 0.23058503 -0.5346335 -0.4963746
## bio7 -0.63911541 -0.3186778 -0.05256181 -0.62275227 -0.3391181 -0.3866963
## bio8 -0.33110726 -0.2212010 0.06324795 -0.37105106 -0.2276076 -0.1386371
## bio9 -0.12411293 -0.6692227 0.60151828 -0.08385793 -0.6317681 -0.6898927
## bio10 -0.36194679 -0.8260330 0.70257598 -0.32811094 -0.8165397 -0.8016505
## bio11 0.09551341 -0.6426766 0.77296096 0.12411688 -0.6183612 -0.5831224
## bio12 0.92211902 0.6216084 -0.20402064 0.91197261 0.6618219 0.6497277
## bio13 1.00000000 0.3184485 0.13929966 0.99092613 0.3638003 0.3707213
## bio14 0.31844853 1.0000000 -0.83281349 0.27096316 0.9913040 0.9791848
## bio15 0.13929966 -0.8328135 1.00000000 0.19326285 -0.8288811 -0.8093866
## bio16 0.99092613 0.2709632 0.19326285 1.00000000 0.3156585 0.3166368
## bio17 0.36380031 0.9913040 -0.82888105 0.31565850 1.0000000 0.9832066
## bio18 0.37072125 0.9791848 -0.80938657 0.31663682 0.9832066 1.0000000
## bio19 0.94856030 0.1661662 0.28320125 0.97382680 0.2114466 0.1873645
## bio19
## bio1 -0.002591276
## bio2 -0.405759856
## bio3 0.296144638
## bio4 -0.589213793
## bio5 -0.327000835
## bio6 0.282322938
## bio7 -0.546340454
## bio8 -0.449569548
## bio9 0.056333201
## bio10 -0.228849357
## bio11 0.183984983
## bio12 0.859723460
## bio13 0.948560299
## bio14 0.166166182
## bio15 0.283201249
## bio16 0.973826798
## bio17 0.211446620
## bio18 0.187364463
## bio19 1.000000000
That’s handy, but maybe not as intuitive as plotting it
raster.cor.plot(spain.worldclim)
## $cor.mds.plot
##
## $cor.heatmap
This contains both a heatmap and an mds plot, let’s look at both!
Okay, now a simpler GLM.
monticola.glm <- enmtools.glm(monticola, spain.worldclim, f = pres ~ poly(bio3, 2) + poly(bio8, 2) + poly(bio11, 2), nback = 400)
interactive.plot.enmtools.model(monticola.glm)
## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO"): Discarded
## ellps WGS 84 in CRS definition: +proj=merc +a=6378137 +b=6378137 +lat_ts=0
## +lon_0=0 +x_0=0 +y_0=0 +k=1 +units=m +nadgrids=@null +wktext +no_defs
## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO"): Discarded
## datum WGS_1984 in CRS definition
## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO"): Discarded
## ellps WGS 84 in CRS definition: +proj=merc +a=6378137 +b=6378137 +lat_ts=0
## +lon_0=0 +x_0=0 +y_0=0 +k=1 +units=m +nadgrids=@null +wktext +no_defs
## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO"): Discarded
## datum WGS_1984 in CRS definition
monticola.glm$response.plots
## $bio3
##
## $bio8
##
## $bio11
Notice that we have flat lines for most of the responses - that’s because we only built a model using variables 1 and 8, but passed it all 19.
Let’s look at what our model looks like in environment space!
visualize.enm(model = monticola.glm, env = spain.worldclim, layers = c("bio8", "bio11"))
## $background.plot
## Warning: Removed 396 rows containing missing values (geom_raster).
##
## $suit.plot
Now let’s try one with some test data for evaluation
monticola.glm <- enmtools.glm(monticola, spain.worldclim, f = pres ~ poly(bio3, 2) + poly(bio8, 2) + poly(bio11, 2), nback = 400, test.prop = 0.3)
monticola.glm
##
##
## Formula: presence ~ poly(bio3, 2) + poly(bio8, 2) + poly(bio11, 2)
## <environment: 0x7fc6cc98cb40>
##
##
## Data table (top ten lines):
##
## | | Longitude| Latitude| bio1| bio2| bio3| bio4| bio5| bio6| bio7| bio8| bio9| bio10| bio11| bio12| bio13| bio14| bio15| bio16| bio17| bio18| bio19| presence|
## |:--|---------:|--------:|----:|----:|----:|----:|----:|----:|----:|----:|----:|-----:|-----:|-----:|-----:|-----:|-----:|-----:|-----:|-----:|-----:|--------:|
## |1 | -6.699292| 43.45172| 121| 88| 41| 4237| 239| 28| 211| 97| 175| 178| 70| 920| 121| 36| 31| 320| 135| 148| 279| 1|
## |2 | -7.188726| 43.43051| 128| 84| 41| 4080| 239| 38| 201| 83| 180| 182| 78| 959| 124| 33| 34| 345| 126| 142| 310| 1|
## |3 | -6.703289| 43.46805| 121| 88| 41| 4237| 239| 28| 211| 97| 175| 178| 70| 920| 121| 36| 31| 320| 135| 148| 279| 1|
## |4 | -6.504768| 43.51037| 134| 85| 42| 3987| 243| 45| 198| 113| 185| 187| 86| 860| 116| 35| 30| 299| 130| 145| 248| 1|
## |6 | -6.728039| 43.55888| 137| 83| 42| 3957| 244| 48| 196| 115| 187| 190| 89| 852| 114| 33| 31| 298| 124| 139| 253| 1|
## |7 | -6.823801| 43.39811| 121| 88| 41| 4237| 239| 28| 211| 97| 175| 178| 70| 920| 121| 36| 31| 320| 135| 148| 279| 1|
## |8 | -6.732560| 43.54083| 137| 83| 42| 3957| 244| 48| 196| 115| 187| 190| 89| 852| 114| 33| 31| 298| 124| 139| 253| 1|
## |10 | -6.541960| 43.45676| 117| 90| 42| 4292| 237| 23| 214| 93| 172| 174| 65| 926| 123| 38| 30| 320| 140| 153| 274| 1|
## |12 | -6.659564| 43.48790| 117| 90| 42| 4292| 237| 23| 214| 93| 172| 174| 65| 926| 123| 38| 30| 320| 140| 153| 274| 1|
## |13 | -6.817537| 43.28005| 105| 93| 41| 4614| 233| 7| 226| 54| 165| 166| 49| 983| 127| 38| 32| 345| 141| 152| 309| 1|
##
##
## Model:
## Call:
## glm(formula = f, family = "binomial", data = analysis.df[, -c(1,
## 2)], weights = weights)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -2.6113 -0.7217 -0.3796 0.6577 2.5628
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -0.1615 0.1067 -1.513 0.130265
## poly(bio3, 2)1 33.4464 3.4375 9.730 < 2e-16 ***
## poly(bio3, 2)2 29.5100 3.8296 7.706 1.30e-14 ***
## poly(bio8, 2)1 -20.3836 3.6370 -5.605 2.09e-08 ***
## poly(bio8, 2)2 7.3182 3.1014 2.360 0.018292 *
## poly(bio11, 2)1 -8.6151 4.1717 -2.065 0.038907 *
## poly(bio11, 2)2 -11.0524 3.2843 -3.365 0.000765 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 815.14 on 693 degrees of freedom
## Residual deviance: 576.05 on 687 degrees of freedom
## AIC: 695.32
##
## Number of Fisher Scoring iterations: 5
##
##
##
## Model fit (training data): class : ModelEvaluation
## n presences : 294
## n absences : 400
## AUC : 0.8546769
## cor : 0.5784453
## max TPR+TNR at : 0.4828914
##
##
## Environment space model fit (training data): class : ModelEvaluation
## n presences : 294
## n absences : 10000
## AUC : 0.6211531
## cor : 0.0859308
## max TPR+TNR at : 0.3840206
##
##
## Proportion of data wittheld for model fitting: 0.3
##
## Model fit (test data): class : ModelEvaluation
## n presences : 127
## n absences : 400
## AUC : 0.8671358
## cor : 0.5276052
## max TPR+TNR at : 0.7213723
##
##
## Environment space model fit (test data): class : ModelEvaluation
## n presences : 127
## n absences : 10000
## AUC : 0.685463
## cor : 0.07766342
## max TPR+TNR at : 0.3840206
##
##
## Suitability:
## class : RasterLayer
## dimensions : 60, 84, 5040 (nrow, ncol, ncell)
## resolution : 0.1666667, 0.1666667 (x, y)
## extent : -10, 4, 35, 45 (xmin, xmax, ymin, ymax)
## crs : +proj=longlat +datum=WGS84 +no_defs
## source : memory
## names : layer
## values : 0.010366, 0.9984437 (min, max)
##
##
##
## Notes:
interactive.plot.enmtools.model(monticola.glm)
## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO"): Discarded
## ellps WGS 84 in CRS definition: +proj=merc +a=6378137 +b=6378137 +lat_ts=0
## +lon_0=0 +x_0=0 +y_0=0 +k=1 +units=m +nadgrids=@null +wktext +no_defs
## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO"): Discarded
## datum WGS_1984 in CRS definition
## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO"): Discarded
## ellps WGS 84 in CRS definition: +proj=merc +a=6378137 +b=6378137 +lat_ts=0
## +lon_0=0 +x_0=0 +y_0=0 +k=1 +units=m +nadgrids=@null +wktext +no_defs
## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO"): Discarded
## datum WGS_1984 in CRS definition
Hey wait! What are these environment space model fit metrics? We’ll get to that in a little while.
monticola.gam <- enmtools.gam(species = monticola, env = spain.worldclim, f = pres ~ bio3 + bio8 + bio11, test.prop = 0.3, nback = 400)
monticola.gam
##
##
## Formula: presence ~ bio3 + bio8 + bio11
## <environment: 0x7fc6de8b7e70>
##
##
## Data table (top ten lines):
##
## | | Longitude| Latitude| bio1| bio2| bio3| bio4| bio5| bio6| bio7| bio8| bio9| bio10| bio11| bio12| bio13| bio14| bio15| bio16| bio17| bio18| bio19| presence|
## |:--|---------:|--------:|----:|----:|----:|----:|----:|----:|----:|----:|----:|-----:|-----:|-----:|-----:|-----:|-----:|-----:|-----:|-----:|-----:|--------:|
## |1 | -6.699292| 43.45172| 121| 88| 41| 4237| 239| 28| 211| 97| 175| 178| 70| 920| 121| 36| 31| 320| 135| 148| 279| 1|
## |3 | -6.703289| 43.46805| 121| 88| 41| 4237| 239| 28| 211| 97| 175| 178| 70| 920| 121| 36| 31| 320| 135| 148| 279| 1|
## |4 | -6.504768| 43.51037| 134| 85| 42| 3987| 243| 45| 198| 113| 185| 187| 86| 860| 116| 35| 30| 299| 130| 145| 248| 1|
## |5 | -6.556106| 43.52801| 134| 85| 42| 3987| 243| 45| 198| 113| 185| 187| 86| 860| 116| 35| 30| 299| 130| 145| 248| 1|
## |6 | -6.728039| 43.55888| 137| 83| 42| 3957| 244| 48| 196| 115| 187| 190| 89| 852| 114| 33| 31| 298| 124| 139| 253| 1|
## |9 | -6.650406| 43.40437| 117| 90| 42| 4292| 237| 23| 214| 93| 172| 174| 65| 926| 123| 38| 30| 320| 140| 153| 274| 1|
## |11 | -6.910225| 43.29794| 117| 91| 41| 4472| 240| 20| 220| 67| 175| 176| 63| 942| 122| 34| 33| 334| 129| 142| 300| 1|
## |12 | -6.659564| 43.48790| 117| 90| 42| 4292| 237| 23| 214| 93| 172| 174| 65| 926| 123| 38| 30| 320| 140| 153| 274| 1|
## |13 | -6.817537| 43.28005| 105| 93| 41| 4614| 233| 7| 226| 54| 165| 166| 49| 983| 127| 38| 32| 345| 141| 152| 309| 1|
## |15 | -6.793121| 43.45933| 121| 88| 41| 4237| 239| 28| 211| 97| 175| 178| 70| 920| 121| 36| 31| 320| 135| 148| 279| 1|
##
##
## Model:
## Family: binomial
## Link function: logit
##
## Formula:
## presence ~ bio3 + bio8 + bio11
##
## Parametric coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -13.038396 1.722752 -7.568 3.78e-14 ***
## bio3 0.389431 0.046443 8.385 < 2e-16 ***
## bio8 -0.030768 0.004762 -6.461 1.04e-10 ***
## bio11 0.004544 0.004534 1.002 0.316
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
##
## R-sq.(adj) = 0.205 Deviance explained = 16.2%
## -REML = 354.25 Scale est. = 1 n = 694
##
##
## Model fit (training data): class : ModelEvaluation
## n presences : 294
## n absences : 400
## AUC : 0.7841667
## cor : 0.4429438
## max TPR+TNR at : 0.5977081
##
##
## Environment space model fit (training data): class : ModelEvaluation
## n presences : 294
## n absences : 10000
## AUC : 0.6670075
## cor : 0.09582997
## max TPR+TNR at : 0.5315484
##
##
## Proportion of data wittheld for model fitting: 0.3
##
## Model fit (test data): class : ModelEvaluation
## n presences : 127
## n absences : 400
## AUC : 0.8291634
## cor : 0.4432638
## max TPR+TNR at : 0.1266629
##
##
## Environment space model fit (test data): class : ModelEvaluation
## n presences : 127
## n absences : 10000
## AUC : 0.7334079
## cor : 0.09095553
## max TPR+TNR at : 0.5316132
##
##
## Suitability:
## class : RasterLayer
## dimensions : 60, 84, 5040 (nrow, ncol, ncell)
## resolution : 0.1666667, 0.1666667 (x, y)
## extent : -10, 4, 35, 45 (xmin, xmax, ymin, ymax)
## crs : +proj=longlat +datum=WGS84 +no_defs
## source : memory
## names : layer
## values : 0.009359476, 0.8559227 (min, max)
##
##
##
## Notes:
## NULL
monticola.gam$response.plots
## $bio3
##
## $bio8
##
## $bio11
visualize.enm(monticola.gam, spain.worldclim, layers = c("bio11", "bio8"))
## $background.plot
## Warning: Removed 396 rows containing missing values (geom_raster).
##
## $suit.plot
If you have maxent installed and setup with dismo, you can even build maxent models.
monticola.mx <- enmtools.maxent(species = monticola, env = spain.worldclim, test.prop = 0.3, nback = 400)
monticola.mx
##
##
## Data table (top ten lines):
##
## | | Longitude| Latitude| presence|
## |:--|---------:|--------:|--------:|
## |1 | -6.699292| 43.45172| 1|
## |2 | -7.188726| 43.43051| 1|
## |3 | -6.703289| 43.46805| 1|
## |4 | -6.504768| 43.51037| 1|
## |7 | -6.823801| 43.39811| 1|
## |9 | -6.650406| 43.40437| 1|
## |11 | -6.910225| 43.29794| 1|
## |12 | -6.659564| 43.48790| 1|
## |13 | -6.817537| 43.28005| 1|
## |14 | -6.831939| 43.38955| 1|
##
##
## Model: Length Class Mode
## 1 MaxEnt S4
##
##
## Model fit (training data): class : ModelEvaluation
## n presences : 294
## n absences : 400
## AUC : 0.9392347
## cor : 0.7782937
## max TPR+TNR at : 0.5464809
##
##
## Environment space model fit (training data): class : ModelEvaluation
## n presences : 294
## n absences : 10000
## AUC : 0.6640782
## cor : 0.1160646
## max TPR+TNR at : 0.5880769
##
##
## Proportion of data wittheld for model fitting: 0.3
##
## Model fit (test data): class : ModelEvaluation
## n presences : 127
## n absences : 400
## AUC : 0.9422244
## cor : 0.7256492
## max TPR+TNR at : 0.4480868
##
##
## Environment space model fit (test data): class : ModelEvaluation
## n presences : 127
## n absences : 10000
## AUC : 0.6878315
## cor : 0.0838308
## max TPR+TNR at : 0.4480868
##
##
## Suitability:
## class : RasterLayer
## dimensions : 60, 84, 5040 (nrow, ncol, ncell)
## resolution : 0.1666667, 0.1666667 (x, y)
## extent : -10, 4, 35, 45 (xmin, xmax, ymin, ymax)
## crs : +proj=longlat +datum=WGS84 +no_defs
## source : memory
## names : layer
## values : 0.004822121, 1 (min, max)
##
##
##
## Notes:
## NULL
monticola.mx$response.plots
## $bio1
##
## $bio2
##
## $bio3
##
## $bio4
##
## $bio5
##
## $bio6
##
## $bio7
##
## $bio8
##
## $bio9
##
## $bio10
##
## $bio11
##
## $bio12
##
## $bio13
##
## $bio14
##
## $bio15
##
## $bio16
##
## $bio17
##
## $bio18
##
## $bio19
visualize.enm(monticola.mx, spain.worldclim, layers = c("bio3", "bio8"))
## $background.plot
## Warning: Removed 396 rows containing missing values (geom_raster).
##
## $suit.plot
Note that what ENMTools is doing for ENM construction is just acting as a front end for dismo and maxent: the interface automates a bunch of stuff, but the underlying math is just dismo and maxent. Stored in the enmtools.model object (monticola.mx), there is a model we can project using dismo’s predict function. We can even use it to bring up the maxent model page in a web browser.
You can even build poisson point process models using ppmlasso. Note: ppmlasso can take a LONG time to fit models with a bunch of predictors, so it’s best to pass a formula here instead of using all 19 variables.
monticola.ppm <- enmtools.ppmlasso(species = monticola, env = spain.worldclim, f = pres ~ bio3 + bio8 + bio11, test.prop = 0.3, nback = 400, gamma = 0)
## Warning in model.matrix.default(mt, mf, contrasts): non-list contrasts argument
## ignored
## Fitting Models: 1 of 200
## Fitting Models: 2 of 200
## Fitting Models: 3 of 200
## Fitting Models: 4 of 200
## Fitting Models: 5 of 200
## Fitting Models: 6 of 200
## Fitting Models: 7 of 200
## Fitting Models: 8 of 200
## Fitting Models: 9 of 200
## Fitting Models: 10 of 200
## Fitting Models: 11 of 200
## Fitting Models: 12 of 200
## Fitting Models: 13 of 200
## Fitting Models: 14 of 200
## Fitting Models: 15 of 200
## Fitting Models: 16 of 200
## Fitting Models: 17 of 200
## Fitting Models: 18 of 200
## Fitting Models: 19 of 200
## Fitting Models: 20 of 200
## Fitting Models: 21 of 200
## Fitting Models: 22 of 200
## Fitting Models: 23 of 200
## Fitting Models: 24 of 200
## Fitting Models: 25 of 200
## Fitting Models: 26 of 200
## Fitting Models: 27 of 200
## Fitting Models: 28 of 200
## Fitting Models: 29 of 200
## Fitting Models: 30 of 200
## Fitting Models: 31 of 200
## Fitting Models: 32 of 200
## Fitting Models: 33 of 200
## Fitting Models: 34 of 200
## Fitting Models: 35 of 200
## Fitting Models: 36 of 200
## Fitting Models: 37 of 200
## Fitting Models: 38 of 200
## Fitting Models: 39 of 200
## Fitting Models: 40 of 200
## Fitting Models: 41 of 200
## Fitting Models: 42 of 200
## Fitting Models: 43 of 200
## Fitting Models: 44 of 200
## Fitting Models: 45 of 200
## Fitting Models: 46 of 200
## Fitting Models: 47 of 200
## Fitting Models: 48 of 200
## Fitting Models: 49 of 200
## Fitting Models: 50 of 200
## Fitting Models: 51 of 200
## Fitting Models: 52 of 200
## Fitting Models: 53 of 200
## Fitting Models: 54 of 200
## Fitting Models: 55 of 200
## Fitting Models: 56 of 200
## Fitting Models: 57 of 200
## Fitting Models: 58 of 200
## Fitting Models: 59 of 200
## Fitting Models: 60 of 200
## Fitting Models: 61 of 200
## Fitting Models: 62 of 200
## Fitting Models: 63 of 200
## Fitting Models: 64 of 200
## Fitting Models: 65 of 200
## Fitting Models: 66 of 200
## Fitting Models: 67 of 200
## Fitting Models: 68 of 200
## Fitting Models: 69 of 200
## Fitting Models: 70 of 200
## Fitting Models: 71 of 200
## Fitting Models: 72 of 200
## Fitting Models: 73 of 200
## Fitting Models: 74 of 200
## Fitting Models: 75 of 200
## Fitting Models: 76 of 200
## Fitting Models: 77 of 200
## Fitting Models: 78 of 200
## Fitting Models: 79 of 200
## Fitting Models: 80 of 200
## Fitting Models: 81 of 200
## Fitting Models: 82 of 200
## Fitting Models: 83 of 200
## Fitting Models: 84 of 200
## Fitting Models: 85 of 200
## Fitting Models: 86 of 200
## Fitting Models: 87 of 200
## Fitting Models: 88 of 200
## Fitting Models: 89 of 200
## Fitting Models: 90 of 200
## Fitting Models: 91 of 200
## Fitting Models: 92 of 200
## Fitting Models: 93 of 200
## Fitting Models: 94 of 200
## Fitting Models: 95 of 200
## Fitting Models: 96 of 200
## Fitting Models: 97 of 200
## Fitting Models: 98 of 200
## Fitting Models: 99 of 200
## Fitting Models: 100 of 200
## Fitting Models: 101 of 200
## Fitting Models: 102 of 200
## Fitting Models: 103 of 200
## Fitting Models: 104 of 200
## Fitting Models: 105 of 200
## Fitting Models: 106 of 200
## Fitting Models: 107 of 200
## Fitting Models: 108 of 200
## Fitting Models: 109 of 200
## Fitting Models: 110 of 200
## Fitting Models: 111 of 200
## Fitting Models: 112 of 200
## Fitting Models: 113 of 200
## Fitting Models: 114 of 200
## Fitting Models: 115 of 200
## Fitting Models: 116 of 200
## Fitting Models: 117 of 200
## Fitting Models: 118 of 200
## Fitting Models: 119 of 200
## Fitting Models: 120 of 200
## Fitting Models: 121 of 200
## Fitting Models: 122 of 200
## Fitting Models: 123 of 200
## Fitting Models: 124 of 200
## Fitting Models: 125 of 200
## Fitting Models: 126 of 200
## Fitting Models: 127 of 200
## Fitting Models: 128 of 200
## Fitting Models: 129 of 200
## Fitting Models: 130 of 200
## Fitting Models: 131 of 200
## Fitting Models: 132 of 200
## Fitting Models: 133 of 200
## Fitting Models: 134 of 200
## Fitting Models: 135 of 200
## Fitting Models: 136 of 200
## Fitting Models: 137 of 200
## Fitting Models: 138 of 200
## Fitting Models: 139 of 200
## Fitting Models: 140 of 200
## Fitting Models: 141 of 200
## Fitting Models: 142 of 200
## Fitting Models: 143 of 200
## Fitting Models: 144 of 200
## Fitting Models: 145 of 200
## Fitting Models: 146 of 200
## Fitting Models: 147 of 200
## Fitting Models: 148 of 200
## Fitting Models: 149 of 200
## Fitting Models: 150 of 200
## Fitting Models: 151 of 200
## Fitting Models: 152 of 200
## Fitting Models: 153 of 200
## Fitting Models: 154 of 200
## Fitting Models: 155 of 200
## Fitting Models: 156 of 200
## Fitting Models: 157 of 200
## Fitting Models: 158 of 200
## Fitting Models: 159 of 200
## Fitting Models: 160 of 200
## Fitting Models: 161 of 200
## Fitting Models: 162 of 200
## Fitting Models: 163 of 200
## Fitting Models: 164 of 200
## Fitting Models: 165 of 200
## Fitting Models: 166 of 200
## Fitting Models: 167 of 200
## Fitting Models: 168 of 200
## Fitting Models: 169 of 200
## Fitting Models: 170 of 200
## Fitting Models: 171 of 200
## Fitting Models: 172 of 200
## Fitting Models: 173 of 200
## Fitting Models: 174 of 200
## Fitting Models: 175 of 200
## Fitting Models: 176 of 200
## Fitting Models: 177 of 200
## Fitting Models: 178 of 200
## Fitting Models: 179 of 200
## Fitting Models: 180 of 200
## Fitting Models: 181 of 200
## Fitting Models: 182 of 200
## Fitting Models: 183 of 200
## Fitting Models: 184 of 200
## Fitting Models: 185 of 200
## Fitting Models: 186 of 200
## Fitting Models: 187 of 200
## Fitting Models: 188 of 200
## Fitting Models: 189 of 200
## Fitting Models: 190 of 200
## Fitting Models: 191 of 200
## Fitting Models: 192 of 200
## Fitting Models: 193 of 200
## Fitting Models: 194 of 200
## Fitting Models: 195 of 200
## Fitting Models: 196 of 200
## Fitting Models: 197 of 200
## Fitting Models: 198 of 200
## Fitting Models: 199 of 200
## Fitting Models: 200 of 200
## Warning in model.matrix.default(mt, mf, contrasts): non-list contrasts argument
## ignored
## Warning in model.matrix.default(mt, mf, contrasts): non-list contrasts argument
## ignored
## Warning in model.matrix.default(mt, mf, contrasts): non-list contrasts argument
## ignored
## Warning in model.matrix.default(mt, mf, contrasts): non-list contrasts argument
## ignored
## Warning in model.matrix.default(mt, mf, contrasts): non-list contrasts argument
## ignored
## Warning in model.matrix.default(mt, mf, contrasts): non-list contrasts argument
## ignored
## Warning in model.matrix.default(mt, mf, contrasts): non-list contrasts argument
## ignored
## Warning in model.matrix.default(mt, mf, contrasts): non-list contrasts argument
## ignored
## Warning in model.matrix.default(mt, mf, contrasts): non-list contrasts argument
## ignored
## Warning in model.matrix.default(mt, mf, contrasts): non-list contrasts argument
## ignored
## Warning in model.matrix.default(mt, mf, contrasts): non-list contrasts argument
## ignored
## Warning in model.matrix.default(mt, mf, contrasts): non-list contrasts argument
## ignored
## Warning in model.matrix.default(mt, mf, contrasts): non-list contrasts argument
## ignored
## Warning in model.matrix.default(mt, mf, contrasts): non-list contrasts argument
## ignored
## Warning in model.matrix.default(mt, mf, contrasts): non-list contrasts argument
## ignored
## Warning in model.matrix.default(mt, mf, contrasts): non-list contrasts argument
## ignored
## Warning in model.matrix.default(mt, mf, contrasts): non-list contrasts argument
## ignored
## Warning in model.matrix.default(mt, mf, contrasts): non-list contrasts argument
## ignored
## Warning in model.matrix.default(mt, mf, contrasts): non-list contrasts argument
## ignored
monticola.ppm
##
##
## Formula: ~bio3 + bio8 + bio11
## <environment: 0x7fc6c52f92d8>
##
##
## Data table (top ten lines):
##
## | | Longitude| Latitude| bio1| bio2| bio3| bio4| bio5| bio6| bio7| bio8| bio9| bio10| bio11| bio12| bio13| bio14| bio15| bio16| bio17| bio18| bio19| presence| wt|
## |:--|---------:|--------:|----:|----:|----:|----:|----:|----:|----:|----:|----:|-----:|-----:|-----:|-----:|-----:|-----:|-----:|-----:|-----:|-----:|--------:|---------:|
## |1 | -6.699292| 43.45172| 121| 88| 41| 4237| 239| 28| 211| 97| 175| 178| 70| 920| 121| 36| 31| 320| 135| 148| 279| 1| 0.0025253|
## |2 | -7.188726| 43.43051| 128| 84| 41| 4080| 239| 38| 201| 83| 180| 182| 78| 959| 124| 33| 34| 345| 126| 142| 310| 1| 0.0069444|
## |4 | -6.504768| 43.51037| 134| 85| 42| 3987| 243| 45| 198| 113| 185| 187| 86| 860| 116| 35| 30| 299| 130| 145| 248| 1| 0.0055556|
## |5 | -6.556106| 43.52801| 134| 85| 42| 3987| 243| 45| 198| 113| 185| 187| 86| 860| 116| 35| 30| 299| 130| 145| 248| 1| 0.0055556|
## |6 | -6.728039| 43.55888| 137| 83| 42| 3957| 244| 48| 196| 115| 187| 190| 89| 852| 114| 33| 31| 298| 124| 139| 253| 1| 0.0092593|
## |7 | -6.823801| 43.39811| 121| 88| 41| 4237| 239| 28| 211| 97| 175| 178| 70| 920| 121| 36| 31| 320| 135| 148| 279| 1| 0.0025253|
## |8 | -6.732560| 43.54083| 137| 83| 42| 3957| 244| 48| 196| 115| 187| 190| 89| 852| 114| 33| 31| 298| 124| 139| 253| 1| 0.0092593|
## |9 | -6.650406| 43.40437| 117| 90| 42| 4292| 237| 23| 214| 93| 172| 174| 65| 926| 123| 38| 30| 320| 140| 153| 274| 1| 0.0055556|
## |10 | -6.541960| 43.45676| 117| 90| 42| 4292| 237| 23| 214| 93| 172| 174| 65| 926| 123| 38| 30| 320| 140| 153| 274| 1| 0.0055556|
## |11 | -6.910225| 43.29794| 117| 91| 41| 4472| 240| 20| 220| 67| 175| 176| 63| 942| 122| 34| 33| 334| 129| 142| 300| 1| 0.0027778|
##
##
## Model: Length Class Mode
## betas 804 -none- numeric
## lambdas 201 -none- numeric
## likelihoods 201 -none- numeric
## pen.likelihoods 201 -none- numeric
## lambda 1 -none- numeric
## beta 4 -none- numeric
## mu 694 -none- numeric
## likelihood 1 -none- numeric
## criterion 1 -none- character
## family 1 -none- character
## gamma 1 -none- numeric
## alpha 1 -none- numeric
## init.coef 1 -none- logical
## criterion.matrix 5 data.frame list
## data 2776 -none- numeric
## pt.interactions 1 -none- logical
## wt 694 -none- numeric
## pres 694 -none- numeric
## x 0 -none- NULL
## y 0 -none- NULL
## r 1 -none- logical
## call 5 -none- call
## formula 2 formula call
## s.means 3 -none- numeric
## s.sds 3 -none- numeric
## cv.group 694 -none- numeric
## n.blocks 1 -none- logical
##
##
## Environment space model fit (training data): class : ModelEvaluation
## n presences : 294
## n absences : 10000
## AUC : 0.6570595
## cor : -0.008231078
## max TPR+TNR at : 24.32807
##
##
## Proportion of data wittheld for model fitting: 0.3
##
## Model fit (test data): class : ModelEvaluation
## n presences : 127
## n absences : 400
## AUC : 0.8027165
## cor : 0.5023459
## max TPR+TNR at : 26.97019
##
##
## Environment space model fit (test data): class : ModelEvaluation
## n presences : 127
## n absences : 10000
## AUC : 0.7078787
## cor : 0.01698327
## max TPR+TNR at : 23.15789
##
##
## Suitability:
## class : RasterLayer
## dimensions : 60, 84, 5040 (nrow, ncol, ncell)
## resolution : 0.1666667, 0.1666667 (x, y)
## extent : -10, 4, 35, 45 (xmin, xmax, ymin, ymax)
## crs : +proj=longlat +datum=WGS84 +no_defs
## source : memory
## names : layer
## values : 0.011549, 2.385701 (min, max)
##
##
##
## Notes:
monticola.ppm$response.plots
## $bio3
##
## $bio8
##
## $bio11
visualize.enm(monticola.ppm, spain.worldclim, layers = c("bio3", "bio8"))
## Warning in model.matrix.default(mt, mf, contrasts): non-list contrasts argument
## ignored
## $background.plot
## Warning: Removed 396 rows containing missing values (geom_raster).
##
## $suit.plot
One final option is to use Blonder et al.’s hypervolume package. At present these are implemented as enmtools.hypervolume objects, rather than as enmtools.model objects. They don’t have as many cool options or visualizations, as many of the options for models don’t make much sense here. It’s worth knowing that you can build them, though, as we will eventually be using them for some of the hypothesis tests.
enmtools.hypervolume(monticola, spain.worldclim[[c(1,8,12,17)]])
##
## Building tree...
## done.
## Ball query...
##
## done.
## Requested probability quantile 0.950000, obtained 0.946393 - setting threshold value 0.000021.
## For a closer match, you can increase num.thresholds in hypervolume_threshold.
## Retaining 203/2037 hypervolume random points for comparison with 5040 test points.
## $hv
## ***** Object of class Hypervolume *****
## Name: Iberolacerta monticola
## Method: Gaussian kernel density estimate
## Number of data points (after weighting): 421
## Dimensionality: 4
## Volume: 15.081957
## Random point density: 135.062046
## Number of random points: 2037
## Random point values:
## min: 0.000
## mean: 0.000
## median: 0.000
## max:0.002
## Parameters:
## kde.bandwidth: 0.2349249 0.2580219 0.2315707 0.1499858
## samples.per.point: 10
## sd.count: 3
## quantile.requested: 0.95
## quantile.requested.type: probability
##
## $suitability
## class : RasterLayer
## dimensions : 60, 84, 5040 (nrow, ncol, ncell)
## resolution : 0.1666667, 0.1666667 (x, y)
## extent : -10, 4, 35, 45 (xmin, xmax, ymin, ymax)
## crs : +proj=longlat +datum=WGS84 +no_defs
## source : memory
## names : bio1
## values : 0, 0.0002667121 (min, max)
##
##
## $species.name
## [1] "Iberolacerta monticola"
##
## $analysis.df
## Longitude Latitude
## 1 -6.699292 43.45172
## 2 -7.188726 43.43051
## 3 -6.703289 43.46805
## 4 -6.504768 43.51037
## 5 -6.556106 43.52801
## 6 -6.728039 43.55888
## 7 -6.823801 43.39811
## 8 -6.732560 43.54083
## 9 -6.650406 43.40437
## 10 -6.541960 43.45676
## 11 -6.910225 43.29794
## 12 -6.659564 43.48790
## 13 -6.817537 43.28005
## 14 -6.831939 43.38955
## 15 -6.793121 43.45933
## 16 -6.499937 43.47067
## 17 -6.619466 43.42050
## 18 -6.724369 43.52515
## 19 -7.657592 40.33878
## 20 -6.657173 43.53675
## 21 -7.627202 40.22996
## 22 -6.661035 43.55847
## 23 -6.773557 43.45048
## 24 -6.673885 43.41994
## 25 -6.660907 43.55608
## 26 -6.933462 43.48396
## 27 -6.687861 43.43833
## 28 -6.951364 43.59288
## 29 -6.688261 43.39621
## 30 -6.704228 43.36865
## 31 -6.619389 43.51894
## 32 -8.150035 43.34445
## 33 -6.618920 43.55074
## 34 -6.770702 43.50165
## 35 -7.754594 40.34196
## 36 -6.770915 43.49353
## 37 -4.896618 43.24043
## 38 -7.654091 40.39647
## 39 -7.922660 43.77674
## 40 -6.746105 43.48776
## 41 -5.975720 43.37362
## 42 -6.958813 43.35228
## 43 -6.604129 43.21463
## 44 -6.674231 43.23169
## 45 -7.049430 42.79107
## 46 -7.154147 42.62575
## 47 -7.128125 42.77214
## 48 -6.231050 43.19318
## 49 -7.613331 40.33810
## 50 -7.613841 40.33794
## 51 -7.587411 40.40216
## 52 -7.056323 43.39149
## 53 -5.621583 43.01901
## 54 -7.982539 43.65810
## 55 -4.869131 43.11211
## 56 -7.810377 43.63642
## 57 -7.824902 43.61708
## 58 -6.904581 43.29112
## 59 -6.827700 42.96513
## 60 -8.040703 43.54231
## 61 -8.008715 43.51093
## 62 -8.067741 43.41750
## 63 -4.980000 43.26000
## 64 -6.931485 43.26939
## 65 -6.763920 43.59477
## 66 -6.927234 43.24961
## 67 -6.180002 43.19246
## 68 -7.407738 43.66977
## 69 -7.407685 43.66973
## 70 -7.407656 43.66979
## 71 -7.407776 43.66982
## 72 -7.161566 42.67164
## 73 -7.740096 43.75711
## 74 -7.480695 40.32510
## 75 -7.729511 40.39481
## 76 -7.613935 40.33855
## 77 -4.853390 43.14163
## 78 -6.264794 43.04339
## 79 -6.276012 43.11746
## 80 -6.964268 43.25487
## 81 -6.745783 43.28473
## 82 -6.898205 43.53368
## 83 -4.599834 43.21588
## 84 -4.851018 43.14197
## 85 -6.139212 43.05136
## 86 -6.142102 42.93155
## 87 -5.171215 43.06957
## 88 -6.036635 43.02531
## 89 -6.651737 43.47214
## 90 -6.690762 43.45555
## 91 -7.942023 43.74054
## 92 -5.402803 42.95040
## 93 -7.473340 43.78935
## 94 -7.411094 43.76633
## 95 -6.575039 42.91070
## 96 -6.137301 43.04235
## 97 -6.238453 43.06750
## 98 -6.139944 43.04458
## 99 -4.890893 43.23017
## 100 -8.180498 43.38845
## 101 -4.931082 43.24231
## 102 -8.063401 43.24602
## 103 -8.053342 43.23926
## 104 -6.834718 42.89244
## 105 -6.580837 42.98877
## 106 -4.878881 43.29832
## 107 -7.599510 40.33032
## 108 -7.597879 40.32729
## 109 -7.571159 40.31599
## 110 -7.563200 40.33297
## 111 -6.053402 43.15157
## 112 -6.050337 43.07146
## 113 -6.132320 43.00186
## 114 -8.144452 43.39255
## 115 -7.589029 40.32022
## 116 -7.618192 40.33799
## 117 -7.787817 40.24217
## 118 -4.686575 43.01294
## 119 -7.588962 40.32031
## 120 -7.588935 40.32903
## 121 -5.048218 43.26821
## 122 -7.668579 40.33195
## 123 -7.582147 40.34275
## 124 -7.595823 40.32304
## 125 -7.589786 40.32417
## 126 -7.641199 40.33058
## 127 -7.590129 40.32201
## 128 -7.150816 43.17367
## 129 -7.612313 40.32155
## 130 -6.876240 42.85325
## 131 -6.760454 42.80724
## 132 -6.854954 42.82537
## 133 -5.132756 43.49572
## 134 -6.329842 43.07315
## 135 -6.230000 43.02000
## 136 -5.972271 43.25346
## 137 -6.225300 43.12128
## 138 -6.207705 43.09089
## 139 -7.787378 40.39362
## 140 -7.647379 40.34105
## 141 -7.570920 40.33198
## 142 -7.572670 40.33843
## 143 -7.605600 40.33840
## 144 -7.340277 42.74906
## 145 -4.812590 43.36524
## 146 -7.146735 42.63721
## 147 -7.136521 42.58219
## 148 -7.597872 40.32673
## 149 -7.587319 40.32879
## 150 -8.110797 43.70124
## 151 -7.612960 40.32116
## 152 -7.614334 40.32299
## 153 -7.610943 40.33872
## 154 -6.114964 43.04907
## 155 -4.941888 43.35310
## 156 -6.971264 43.50353
## 157 -7.596052 40.32335
## 158 -7.645674 40.36543
## 159 -7.621731 40.34170
## 160 -7.642539 40.36317
## 161 -6.803279 43.30922
## 162 -6.968332 43.34016
## 163 -6.961694 43.34072
## 164 -6.644926 43.04650
## 165 -7.482150 40.38189
## 166 -4.780000 43.02000
## 167 -4.780000 43.11000
## 168 -4.780000 43.20000
## 169 -4.790000 43.29000
## 170 -4.910000 43.11000
## 171 -4.910000 43.20000
## 172 -5.030000 43.11000
## 173 -5.030000 43.20000
## 174 -5.150000 43.11000
## 175 -5.270000 43.11000
## 176 -5.390000 43.01000
## 177 -5.400000 43.10000
## 178 -5.520000 43.01000
## 179 -5.520000 43.10000
## 180 -5.760000 43.00000
## 181 -5.880000 43.00000
## 182 -6.110000 43.09000
## 183 -6.120000 43.00000
## 184 -6.240000 43.00000
## 185 -6.240000 43.09000
## 186 -6.360000 43.10000
## 187 -6.520000 42.20000
## 188 -6.610000 42.92000
## 189 -6.640000 42.20000
## 190 -6.760000 42.21000
## 191 -6.760000 42.30000
## 192 -6.850000 43.02000
## 193 -6.850000 43.20000
## 194 -6.860000 42.75000
## 195 -6.860000 42.84000
## 196 -6.860000 42.93000
## 197 -6.870000 42.39000
## 198 -6.880000 42.12000
## 199 -6.970000 43.20000
## 200 -6.980000 42.75000
## 201 -6.980000 43.02000
## 202 -6.990000 42.57000
## 203 -6.990000 42.66000
## 204 -7.010000 41.94000
## 205 -7.080000 43.47000
## 206 -7.090000 43.20000
## 207 -7.090000 43.29000
## 208 -7.090000 43.38000
## 209 -7.110000 42.57000
## 210 -7.210000 43.38000
## 211 -7.240000 42.30000
## 212 -7.700000 43.48000
## 213 -7.700000 43.66000
## 214 -7.700000 43.75000
## 215 -7.830000 43.48000
## 216 -7.950000 43.30000
## 217 -7.950000 43.39000
## 218 -7.960000 43.03000
## 219 -7.820000 43.57000
## 220 -7.820000 43.66000
## 221 -7.820000 43.75000
## 222 -7.940000 43.75000
## 223 -7.950000 43.48000
## 224 -7.950000 43.57000
## 225 -7.950000 43.66000
## 226 -8.070000 43.39000
## 227 -8.070000 43.48000
## 228 -8.070000 43.66000
## 229 -8.070000 43.75000
## 230 -8.080000 43.03000
## 231 -8.080000 43.21000
## 232 -8.080000 43.30000
## 233 -8.200000 43.31000
## 234 -8.200000 43.40000
## 235 -8.320000 43.31000
## 236 -9.060000 42.86000
## 237 -5.150000 43.02000
## 238 -5.270000 43.02000
## 239 -5.390000 42.92000
## 240 -5.510000 42.83000
## 241 -5.630000 42.92000
## 242 -6.390000 42.29000
## 243 -6.610000 42.83000
## 244 -6.630000 42.29000
## 245 -6.740000 42.75000
## 246 -6.740000 42.84000
## 247 -6.980000 42.84000
## 248 -7.090000 43.11000
## 249 -7.100000 42.84000
## 250 -7.100000 42.93000
## 251 -7.100000 43.02000
## 252 -7.110000 42.66000
## 253 -7.210000 43.29000
## 254 -7.330000 43.39000
## 255 -7.330000 43.48000
## 256 -7.350000 42.67000
## 257 -7.450000 43.48000
## 258 -7.570000 43.66000
## 259 -7.580000 43.48000
## 260 -7.240000 42.21000
## 261 -7.360000 42.22000
## 262 -7.490000 42.22000
## 263 -4.910000 43.29000
## 264 -5.030000 43.29000
## 265 -5.160000 43.47000
## 266 -5.280000 43.20000
## 267 -5.290000 43.47000
## 268 -5.400000 43.19000
## 269 -5.520000 43.19000
## 270 -5.890000 43.09000
## 271 -5.890000 43.18000
## 272 -6.110000 43.18000
## 273 -6.480000 43.10000
## 274 -6.600000 43.28000
## 275 -6.610000 43.01000
## 276 -6.710000 43.47000
## 277 -6.720000 43.20000
## 278 -6.720000 43.29000
## 279 -6.730000 43.02000
## 280 -6.730000 43.11000
## 281 -6.840000 43.38000
## 282 -6.840000 43.47000
## 283 -6.960000 43.38000
## 284 -6.960000 43.47000
## 285 -6.960000 43.56000
## 286 -6.970000 43.29000
## 287 -6.760000 42.12000
## 288 -6.121960 43.05522
## 289 -6.214485 43.07152
## 290 -5.005530 43.26170
## 291 -5.820000 43.14000
## 292 -7.613616 40.32185
## 293 -4.910000 43.15000
## 294 -7.363881 43.42967
## 295 -6.430000 42.88000
## 296 -6.190000 42.85000
## 297 -6.550000 43.18000
## 298 -6.660000 42.11000
## 299 -7.616667 40.31667
## 300 -7.638880 40.36111
## 301 -6.440000 42.26000
## 302 -6.730000 42.82000
## 303 -6.620000 42.77000
## 304 -6.580000 42.27000
## 305 -5.530000 42.89000
## 306 -0.020000 42.59000
## 307 -2.450000 42.46000
## 308 -7.630000 40.33000
## 309 -8.080000 43.03000
## 310 -8.320000 43.31000
## 311 -9.060000 42.86000
## 312 -8.200000 43.31000
## 313 -8.080000 43.21000
## 314 -7.960000 43.03000
## 315 -7.950000 43.30000
## 316 -8.080000 43.30000
## 317 -7.090000 43.29000
## 318 -6.970000 43.20000
## 319 -7.950000 43.39000
## 320 -8.070000 43.39000
## 321 -8.200000 43.40000
## 322 -7.090000 43.20000
## 323 -7.090000 43.11000
## 324 -6.970000 43.29000
## 325 -7.240000 42.21000
## 326 -6.630000 42.29000
## 327 -7.210000 43.29000
## 328 -7.490000 42.22000
## 329 -7.350000 42.67000
## 330 -6.870000 42.39000
## 331 -7.240000 42.30000
## 332 -7.360000 42.22000
## 333 -6.760000 42.30000
## 334 -6.760000 42.21000
## 335 -6.760000 42.12000
## 336 -7.010000 41.94000
## 337 -6.640000 42.20000
## 338 -6.880000 42.12000
## 339 -5.630000 42.92000
## 340 -4.830000 40.50000
## 341 -5.510000 42.83000
## 342 -5.510000 42.92000
## 343 -4.710000 40.32000
## 344 -5.650000 40.31000
## 345 -6.340000 43.46000
## 346 -5.520000 43.01000
## 347 -5.520000 43.19000
## 348 -5.970000 42.99000
## 349 -5.770000 40.30000
## 350 -5.520000 43.10000
## 351 -5.180000 40.41000
## 352 -5.970000 43.09000
## 353 -4.820000 40.32000
## 354 -6.610000 42.92000
## 355 -7.200000 43.56000
## 356 -6.480000 43.10000
## 357 -6.720000 43.20000
## 358 -6.360000 43.01000
## 359 -6.720000 43.29000
## 360 -6.610000 42.83000
## 361 -6.610000 43.01000
## 362 -6.370000 42.92000
## 363 -6.710000 43.47000
## 364 -7.580000 43.57000
## 365 -6.840000 43.38000
## 366 -6.850000 43.20000
## 367 -6.960000 43.56000
## 368 -7.330000 43.39000
## 369 -7.210000 43.38000
## 370 -7.080000 43.47000
## 371 -6.730000 43.11000
## 372 -6.960000 43.47000
## 373 -6.520000 42.20000
## 374 -6.840000 43.47000
## 375 -6.110000 40.57000
## 376 -7.090000 43.38000
## 377 -6.390000 42.29000
## 378 -5.990000 43.27000
## 379 -6.360000 43.10000
## 380 -6.600000 43.28000
## 381 -6.350000 43.19000
## 382 -6.960000 43.38000
## 383 -6.010000 43.13000
## 384 -6.020000 42.99000
## 385 -4.910000 43.11000
## 386 -5.150000 43.11000
## 387 -4.910000 43.20000
## 388 0.130000 42.64000
## 389 -5.390000 42.92000
## 390 -5.280000 43.20000
## 391 -3.890000 40.87000
## 392 -5.030000 43.11000
## 393 -5.030000 43.20000
## 394 -4.780000 43.11000
## 395 -5.400000 43.10000
## 396 -4.910000 43.29000
## 397 -5.390000 43.01000
## 398 -4.130000 40.78000
## 399 -3.890000 40.96000
## 400 -3.770000 40.78000
## 401 -5.030000 43.29000
## 402 -4.790000 43.29000
## 403 -5.270000 43.11000
## 404 -4.660000 43.21000
## 405 -3.890000 40.69000
## 406 -4.590000 40.32000
## 407 -4.780000 43.20000
## 408 -5.160000 43.47000
## 409 -3.770000 41.05000
## 410 -5.290000 43.47000
## 411 -5.400000 43.19000
## 412 -4.010000 40.69000
## 413 -4.010000 40.78000
## 414 -3.770000 40.96000
## 415 -4.650000 42.94000
## 416 -5.280000 43.38000
## 417 -3.890000 40.78000
## 418 -4.120000 40.60000
## 419 -5.160000 43.38000
## 420 -3.650000 41.05000
## 421 -4.010000 40.87000
##
## attr(,"class")
## [1] "enmtools.hypervolume"
Often we have applications where we want to take a model made in one time and place and project it to the distribution of environments in another time and place. This is relatively easy to do in ENMTools! For instance we can extrapolate our model across the entire world.
gam.pred <- predict(monticola.gam, all.worldclim)
gam.pred$suitability.plot
## Warning: Raster pixels are placed at uneven vertical intervals and will be
## shifted. Consider using geom_tile() instead.
It’s a bit crazy to think that these lizards are going to live at the north pole, so let’s look at the built-in plots that visualize the environments where the model was trained vs. where it’s being projected.
gam.pred$threespace.plot
Wanna know where your model’s predictions are being reigned in by clamping? ENMTools has got you covered.
gam.pred$clamp.plot
## Warning: Raster pixels are placed at uneven vertical intervals and will be
## shifted. Consider using geom_tile() instead.
We can look at the suitability projected in the future as well. Note: this involves downloading new layers so this bit of code may take a while to run.
spain.future <- raster::getData('CMIP5', var = "bio", res = 10, rcp = 85, model = 'AC', year = 70)
## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO"): Discarded datum Unknown based on WGS84 ellipsoid in CRS definition,
## but +towgs84= values preserved
## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO"): Discarded datum Unknown based on WGS84 ellipsoid in CRS definition,
## but +towgs84= values preserved
## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO"): Discarded datum Unknown based on WGS84 ellipsoid in CRS definition,
## but +towgs84= values preserved
## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO"): Discarded datum Unknown based on WGS84 ellipsoid in CRS definition,
## but +towgs84= values preserved
## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO"): Discarded datum Unknown based on WGS84 ellipsoid in CRS definition,
## but +towgs84= values preserved
## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO"): Discarded datum Unknown based on WGS84 ellipsoid in CRS definition,
## but +towgs84= values preserved
## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO"): Discarded datum Unknown based on WGS84 ellipsoid in CRS definition,
## but +towgs84= values preserved
## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO"): Discarded datum Unknown based on WGS84 ellipsoid in CRS definition,
## but +towgs84= values preserved
## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO"): Discarded datum Unknown based on WGS84 ellipsoid in CRS definition,
## but +towgs84= values preserved
## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO"): Discarded datum Unknown based on WGS84 ellipsoid in CRS definition,
## but +towgs84= values preserved
## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO"): Discarded datum Unknown based on WGS84 ellipsoid in CRS definition,
## but +towgs84= values preserved
## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO"): Discarded datum Unknown based on WGS84 ellipsoid in CRS definition,
## but +towgs84= values preserved
## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO"): Discarded datum Unknown based on WGS84 ellipsoid in CRS definition,
## but +towgs84= values preserved
## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO"): Discarded datum Unknown based on WGS84 ellipsoid in CRS definition,
## but +towgs84= values preserved
## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO"): Discarded datum Unknown based on WGS84 ellipsoid in CRS definition,
## but +towgs84= values preserved
## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO"): Discarded datum Unknown based on WGS84 ellipsoid in CRS definition,
## but +towgs84= values preserved
## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO"): Discarded datum Unknown based on WGS84 ellipsoid in CRS definition,
## but +towgs84= values preserved
## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO"): Discarded datum Unknown based on WGS84 ellipsoid in CRS definition,
## but +towgs84= values preserved
## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO"): Discarded datum Unknown based on WGS84 ellipsoid in CRS definition,
## but +towgs84= values preserved
## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO"): Discarded datum Unknown based on WGS84 ellipsoid in CRS definition,
## but +towgs84= values preserved
spain.extent <- extent(-10, 4, 35, 45)
spain.future <- crop(x = spain.future, y = spain.extent)
names(spain.future) <- names(spain.worldclim)
future.pred <- predict(monticola.gam, spain.future)
future.pred$suitability.plot
future.pred$threespace.plot
future.pred$clamp.plot
Note that you can call predict on the model objects stored inside ENMTools models as well, but you often need to change the order. For instance to predict using a maxent model, you pass in the environmental rasters first and the model second. This is just a dismo thing to be aware of, there’s no particular logic to either ordering that you need to understand.
plot(predict(all.worldclim, monticola.gam$model))